home *** CD-ROM | disk | FTP | other *** search
/ Aminet 1 (Walnut Creek) / Aminet - June 1993 [Walnut Creek].iso / usenet / sources / volume2 / aplictns / matlab / src.7 < prev   
Internet Message Format  |  1988-11-02  |  15KB
  1. Path: xanth!nic.MR.NET!hal!cwjcc!mailrus!bbn!ulowell!page
  2. From: page@swan.ulowell.edu (Bob Page)
  3. Newsgroups: comp.sources.amiga
  4. Subject: v02i047:  matlab - matrix laboratory, Part07/11
  5. Message-ID: <10022@swan.ulowell.edu>
  6. Date: 2 Nov 88 21:46:08 GMT
  7. Organization: University of Lowell, Computer Science Dept.
  8. Lines: 388
  9. Approved: page@swan.ulowell.edu
  10.  
  11. Submitted-by: strovink%galaxy-43@afit-ab.arpa (Mark A. Strovink)
  12. Posting-number: Volume 2, Issue 47
  13. Archive-name: applications/matlab/src.7
  14.  
  15. #    This is a shell archive.
  16. #    Remove everything above and including the cut line.
  17. #    Then run the rest of the file through sh.
  18. #----cut here-----cut here-----cut here-----cut here----#
  19. #!/bin/sh
  20. # shar:    Shell Archiver
  21. #    Run the following text with /bin/sh to create:
  22. #    src-7
  23. # This archive created: Wed Nov  2 16:20:51 1988
  24. cat << \SHAR_EOF > src-7
  25.       CR = FLOP((ARS*BRS + AIS*BIS)/D)    
  26.       CI = (AIS*BRS - ARS*BIS)/D          
  27.       IF (CI .NE. 0.0D0) CI = FLOP(CI)    
  28.       RETURN           
  29.       END              
  30.       SUBROUTINE WSIGN(XR,XI,YR,YI,ZR,ZI)                    
  31.       DOUBLE PRECISION XR,XI,YR,YI,ZR,ZI,PYTHAG,T            
  32. C     IF Y .NE. 0, Z = X*Y/ABS(Y)         
  33. C     IF Y .EQ. 0, Z = X                  
  34.       T = PYTHAG(YR,YI)                   
  35.       ZR = XR          
  36.       ZI = XI          
  37.       IF (T .NE. 0.0D0) CALL WMUL(YR/T,YI/T,ZR,ZI,ZR,ZI)     
  38.       RETURN           
  39.       END              
  40.       SUBROUTINE WSQRT(XR,XI,YR,YI)       
  41.       DOUBLE PRECISION XR,XI,YR,YI,S,TR,TI,PYTHAG,FLOP       
  42. C     Y = SQRT(X) WITH YR .GE. 0.0 AND SIGN(YI) .EQ. SIGN(XI)                   
  43. C                      
  44.       TR = XR          
  45.       TI = XI          
  46.       S = DSQRT(0.5D0*(PYTHAG(TR,TI) + DABS(TR)))            
  47.       IF (TR .GE. 0.0D0) YR = FLOP(S)     
  48.       IF (TI .LT. 0.0D0) S = -S           
  49.       IF (TR .LE. 0.0D0) YI = FLOP(S)     
  50.       IF (TR .LT. 0.0D0) YR = FLOP(0.5D0*(TI/YI))            
  51.       IF (TR .GT. 0.0D0) YI = FLOP(0.5D0*(TI/YR))            
  52.       RETURN           
  53.       END              
  54.       SUBROUTINE WLOG(XR,XI,YR,YI)        
  55.       DOUBLE PRECISION XR,XI,YR,YI,T,R,PYTHAG                
  56. C     Y = LOG(X)       
  57.       R = PYTHAG(XR,XI)                   
  58.       IF (R .EQ. 0.0D0) CALL ERROR(32)    
  59.       IF (R .EQ. 0.0D0) RETURN            
  60.       T = DATAN2(XI,XR)                   
  61.       IF (XI.EQ.0.0D0 .AND. XR.LT.0.0D0) T = DABS(T)         
  62.       YR = DLOG(R)     
  63.       YI = T           
  64.       RETURN           
  65.       END              
  66.       SUBROUTINE WATAN(XR,XI,YR,YI)       
  67. C     Y = ATAN(X) = (I/2)*LOG((I+X)/(I-X))                   
  68.       DOUBLE PRECISION XR,XI,YR,YI,TR,TI  
  69.       IF (XI .NE. 0.0D0) GO TO 10         
  70.          YR = DATAN2(XR,1.0D0)            
  71.          YI = 0.0D0    
  72.          RETURN        
  73.    10 IF (XR.NE.0.0D0 .OR. DABS(XI).NE.1.0D0) GO TO 20       
  74.          CALL ERROR(32)                   
  75.          RETURN        
  76.    20 CALL WDIV(XR,1.0D0+XI,-XR,1.0D0-XI,TR,TI)              
  77.       CALL WLOG(TR,TI,TR,TI)              
  78.       YR = -TI/2.0D0   
  79.       YI = TR/2.0D0    
  80.       RETURN           
  81.       END              
  82.       DOUBLE PRECISION FUNCTION WNRM2(N,XR,XI,INCX)          
  83.       DOUBLE PRECISION XR(1),XI(1),PYTHAG,S                  
  84. C     NORM2(X)         
  85.       S = 0.0D0        
  86.       IF (N .LE. 0) GO TO 20              
  87.       IX = 1           
  88.       DO 10 I = 1, N   
  89.          S = PYTHAG(S,XR(IX))             
  90.          S = PYTHAG(S,XI(IX))             
  91.          IX = IX + INCX                   
  92.    10 CONTINUE         
  93.    20 WNRM2 = S        
  94.       RETURN           
  95.       END              
  96.       DOUBLE PRECISION FUNCTION WASUM(N,XR,XI,INCX)          
  97.       DOUBLE PRECISION XR(1),XI(1),S,FLOP                    
  98. C     NORM1(X)         
  99.       S = 0.0D0        
  100.       IF (N .LE. 0) GO TO 20              
  101.       IX = 1           
  102.       DO 10 I = 1, N   
  103.          S = FLOP(S + DABS(XR(IX)) + DABS(XI(IX)))           
  104.          IX = IX + INCX                   
  105.    10 CONTINUE         
  106.    20 WASUM = S        
  107.       RETURN           
  108.       END              
  109.       INTEGER FUNCTION IWAMAX(N,XR,XI,INCX)                  
  110.       DOUBLE PRECISION XR(1),XI(1),S,P    
  111. C     INDEX OF NORMINF(X)                 
  112.       K = 0            
  113.       IF (N .LE. 0) GO TO 20              
  114.       K = 1            
  115.       S = 0.0D0        
  116.       IX = 1           
  117.       DO 10 I = 1, N   
  118.          P = DABS(XR(IX)) + DABS(XI(IX))  
  119.          IF (P .GT. S) K = I              
  120.          IF (P .GT. S) S = P              
  121.          IX = IX + INCX                   
  122.    10 CONTINUE         
  123.    20 IWAMAX = K       
  124.       RETURN           
  125.       END              
  126.       SUBROUTINE WRSCAL(N,S,XR,XI,INCX)   
  127.       DOUBLE PRECISION S,XR(1),XI(1),FLOP                    
  128.       IF (N .LE. 0) RETURN                
  129.       IX = 1           
  130.       DO 10 I = 1, N   
  131.          XR(IX) = FLOP(S*XR(IX))          
  132.          IF (XI(IX) .NE. 0.0D0) XI(IX) = FLOP(S*XI(IX))      
  133.          IX = IX + INCX                   
  134.    10 CONTINUE         
  135.       RETURN           
  136.       END              
  137.       SUBROUTINE WSCAL(N,SR,SI,XR,XI,INCX)                   
  138.       DOUBLE PRECISION SR,SI,XR(1),XI(1)  
  139.       IF (N .LE. 0) RETURN                
  140.       IX = 1           
  141.       DO 10 I = 1, N   
  142.          CALL WMUL(SR,SI,XR(IX),XI(IX),XR(IX),XI(IX))        
  143.          IX = IX + INCX                   
  144.    10 CONTINUE         
  145.       RETURN           
  146.       END              
  147.       SUBROUTINE WAXPY(N,SR,SI,XR,XI,INCX,YR,YI,INCY)        
  148.       DOUBLE PRECISION SR,SI,XR(1),XI(1),YR(1),YI(1),FLOP    
  149.       IF (N .LE. 0) RETURN                
  150.       IF (SR .EQ. 0.0D0 .AND. SI .EQ. 0.0D0) RETURN          
  151.       IX = 1           
  152.       IY = 1           
  153.       IF (INCX.LT.0) IX = (-N+1)*INCX + 1                    
  154.       IF (INCY.LT.0) IY = (-N+1)*INCY + 1                    
  155.       DO 10 I = 1, N   
  156.          YR(IY) = FLOP(YR(IY) + SR*XR(IX) - SI*XI(IX))       
  157.          YI(IY) = YI(IY) + SR*XI(IX) + SI*XR(IX)             
  158.          IF (YI(IY) .NE. 0.0D0) YI(IY) = FLOP(YI(IY))        
  159.          IX = IX + INCX                   
  160.          IY = IY + INCY                   
  161.    10 CONTINUE         
  162.       RETURN           
  163.       END              
  164.       DOUBLE PRECISION FUNCTION WDOTUR(N,XR,XI,INCX,YR,YI,INCY)                 
  165.       DOUBLE PRECISION XR(1),XI(1),YR(1),YI(1),S,FLOP        
  166.       S = 0.0D0        
  167.       IF (N .LE. 0) GO TO 20              
  168.       IX = 1           
  169.       IY = 1           
  170.       IF (INCX.LT.0) IX = (-N+1)*INCX + 1                    
  171.       IF (INCY.LT.0) IY = (-N+1)*INCY + 1                    
  172.       DO 10 I = 1, N   
  173.          S = FLOP(S + XR(IX)*YR(IY) - XI(IX)*YI(IY))         
  174.          IX = IX + INCX                   
  175.          IY = IY + INCY                   
  176.    10 CONTINUE         
  177.    20 WDOTUR = S       
  178.       RETURN           
  179.       END              
  180.       DOUBLE PRECISION FUNCTION WDOTUI(N,XR,XI,INCX,YR,YI,INCY)                 
  181.       DOUBLE PRECISION XR(1),XI(1),YR(1),YI(1),S,FLOP        
  182.       S = 0.0D0        
  183.       IF (N .LE. 0) GO TO 20              
  184.       IX = 1           
  185.       IY = 1           
  186.       IF (INCX.LT.0) IX = (-N+1)*INCX + 1                    
  187.       IF (INCY.LT.0) IY = (-N+1)*INCY + 1                    
  188.       DO 10 I = 1, N   
  189.          S = S + XR(IX)*YI(IY) + XI(IX)*YR(IY)               
  190.          IF (S .NE. 0.0D0) S = FLOP(S)    
  191.          IX = IX + INCX                   
  192.          IY = IY + INCY                   
  193.    10 CONTINUE         
  194.    20 WDOTUI = S       
  195.       RETURN           
  196.       END              
  197.       DOUBLE PRECISION FUNCTION WDOTCR(N,XR,XI,INCX,YR,YI,INCY)                 
  198.       DOUBLE PRECISION XR(1),XI(1),YR(1),YI(1),S,FLOP        
  199.       S = 0.0D0        
  200.       IF (N .LE. 0) GO TO 20              
  201.       IX = 1           
  202.       IY = 1           
  203.       IF (INCX.LT.0) IX = (-N+1)*INCX + 1                    
  204.       IF (INCY.LT.0) IY = (-N+1)*INCY + 1                    
  205.       DO 10 I = 1, N   
  206.          S = FLOP(S + XR(IX)*YR(IY) + XI(IX)*YI(IY))         
  207.          IX = IX + INCX                   
  208.          IY = IY + INCY                   
  209.    10 CONTINUE         
  210.    20 WDOTCR = S       
  211.       RETURN           
  212.       END              
  213.       DOUBLE PRECISION FUNCTION WDOTCI(N,XR,XI,INCX,YR,YI,INCY)                 
  214.       DOUBLE PRECISION XR(1),XI(1),YR(1),YI(1),S,FLOP        
  215.       S = 0.0D0        
  216.       IF (N .LE. 0) GO TO 20              
  217.       IX = 1           
  218.       IY = 1           
  219.       IF (INCX.LT.0) IX = (-N+1)*INCX + 1                    
  220.       IF (INCY.LT.0) IY = (-N+1)*INCY + 1                    
  221.       DO 10 I = 1, N   
  222.          S = S + XR(IX)*YI(IY) - XI(IX)*YR(IY)               
  223.          IF (S .NE. 0.0D0) S = FLOP(S)    
  224.          IX = IX + INCX                   
  225.          IY = IY + INCY                   
  226.    10 CONTINUE         
  227.    20 WDOTCI = S       
  228.       RETURN           
  229.       END              
  230.       SUBROUTINE WCOPY(N,XR,XI,INCX,YR,YI,INCY)              
  231.       DOUBLE PRECISION XR(1),XI(1),YR(1),YI(1)               
  232.       IF (N .LE. 0) RETURN                
  233.       IX = 1           
  234.       IY = 1           
  235.       IF (INCX.LT.0) IX = (-N+1)*INCX + 1                    
  236.       IF (INCY.LT.0) IY = (-N+1)*INCY + 1                    
  237.       DO 10 I = 1, N   
  238.          YR(IY) = XR(IX)                  
  239.          YI(IY) = XI(IX)                  
  240.          IX = IX + INCX                   
  241.          IY = IY + INCY                   
  242.    10 CONTINUE         
  243.       RETURN           
  244.       END              
  245.       SUBROUTINE WSET(N,XR,XI,YR,YI,INCY)                    
  246.       INTEGER N,INCY   
  247.       DOUBLE PRECISION XR,XI,YR(1),YI(1)  
  248.       IY = 1           
  249.       IF (N .LE. 0 ) RETURN               
  250.       DO 10 I = 1,N    
  251.          YR(IY) = XR   
  252.          YI(IY) = XI   
  253.          IY = IY + INCY                   
  254.    10 CONTINUE         
  255.       RETURN           
  256.       END              
  257.       SUBROUTINE WSWAP(N,XR,XI,INCX,YR,YI,INCY)              
  258.       DOUBLE PRECISION XR(1),XI(1),YR(1),YI(1),T             
  259.       IF (N .LE. 0) RETURN                
  260.       IX = 1           
  261.       IY = 1           
  262.       IF (INCX.LT.0) IX = (-N+1)*INCX + 1                    
  263.       IF (INCY.LT.0) IY = (-N+1)*INCY + 1                    
  264.       DO 10 I = 1, N   
  265.          T = XR(IX)    
  266.          XR(IX) = YR(IY)                  
  267.          YR(IY) = T    
  268.          T = XI(IX)    
  269.          XI(IX) = YI(IY)                  
  270.          YI(IY) = T    
  271.          IX = IX + INCX                   
  272.          IY = IY + INCY                   
  273.    10 CONTINUE         
  274.       RETURN           
  275.       END              
  276.       SUBROUTINE RSET(N,DX,DY,INCY)       
  277. C                      
  278. C     COPIES A SCALAR, X, TO A SCALAR, Y.                    
  279.       DOUBLE PRECISION DX,DY(1)           
  280. C                      
  281.       IF (N.LE.0) RETURN                  
  282.       IY = 1           
  283.       IF (INCY.LT.0) IY = (-N+1)*INCY + 1                    
  284.       DO 10 I = 1,N    
  285.         DY(IY) = DX    
  286.         IY = IY + INCY                    
  287.    10 CONTINUE         
  288.       RETURN           
  289.       END              
  290.       SUBROUTINE RSWAP(N,X,INCX,Y,INCY)   
  291.       DOUBLE PRECISION X(1),Y(1),T        
  292.       IF (N .LE. 0) RETURN                
  293.       IX = 1           
  294.       IY = 1           
  295.       IF (INCX.LT.0) IX = (-N+1)*INCX+1   
  296.       IF (INCY.LT.0) IY = (-N+1)*INCY+1   
  297.       DO 10 I = 1, N   
  298.          T = X(IX)     
  299.          X(IX) = Y(IY)                    
  300.          Y(IY) = T     
  301.          IX = IX + INCX                   
  302.          IY = IY + INCY                   
  303.    10 CONTINUE         
  304.       RETURN           
  305.       END              
  306.       SUBROUTINE RROT(N,DX,INCX,DY,INCY,C,S)                 
  307. C                      
  308. C     APPLIES A PLANE ROTATION.           
  309.       DOUBLE PRECISION DX(1),DY(1),DTEMP,C,S,FLOP            
  310.       INTEGER I,INCX,INCY,IX,IY,N         
  311. C                      
  312.       IF (N.LE.0) RETURN                  
  313.       IX = 1           
  314.       IY = 1           
  315.       IF (INCX.LT.0) IX = (-N+1)*INCX + 1                    
  316.       IF (INCY.LT.0) IY = (-N+1)*INCY + 1                    
  317.       DO 10 I = 1,N    
  318.         DTEMP = FLOP(C*DX(IX) + S*DY(IY))                    
  319.         DY(IY) = FLOP(C*DY(IY) - S*DX(IX))                   
  320.         DX(IX) = DTEMP                    
  321.         IX = IX + INCX                    
  322.         IY = IY + INCY                    
  323.    10 CONTINUE         
  324.       RETURN           
  325.       END              
  326.       SUBROUTINE RROTG(DA,DB,C,S)         
  327. C                      
  328. C     CONSTRUCT GIVENS PLANE ROTATION.    
  329. C                      
  330.       DOUBLE PRECISION DA,DB,C,S,RHO,PYTHAG,FLOP,R,Z         
  331. C                      
  332.       RHO = DB         
  333.       IF ( DABS(DA) .GT. DABS(DB) ) RHO = DA                 
  334.       C = 1.0D0        
  335.       S = 0.0D0        
  336.       Z = 1.0D0        
  337.       R = FLOP(DSIGN(PYTHAG(DA,DB),RHO))  
  338.       IF (R .NE. 0.0D0) C = FLOP(DA/R)    
  339.       IF (R .NE. 0.0D0) S = FLOP(DB/R)    
  340.       IF ( DABS(DA) .GT. DABS(DB) ) Z = S                    
  341.       IF ( DABS(DB) .GE. DABS(DA) .AND. C .NE. 0.0D0 ) Z = FLOP(1.0D0/C)        
  342.       DA = R           
  343.       DB = Z           
  344.       RETURN           
  345.       END              
  346.       LOGICAL FUNCTION EQID(X,Y)          
  347. C     CHECK FOR EQUALITY OF TWO NAMES     
  348.       INTEGER X(4),Y(4)                   
  349.       EQID = .TRUE.    
  350.       DO 10 I = 1, 4   
  351.    10 EQID = EQID .AND. (X(I).EQ.Y(I))    
  352.       RETURN           
  353.       END              
  354.       SUBROUTINE PUTID(X,Y)               
  355. C     STORE A NAME     
  356.       INTEGER X(4),Y(4)                   
  357.       DO 10 I = 1, 4   
  358.    10 X(I) = Y(I)      
  359.       RETURN           
  360.       END              
  361.       DOUBLE PRECISION FUNCTION ROUND(X)  
  362.       DOUBLE PRECISION X,Y,Z,E,H          
  363.       DATA H/1.0D9/    
  364.       Z = DABS(X)      
  365.       Y = Z + 1.0D0    
  366.       IF (Y .EQ. Z) GO TO 40              
  367.       Y = 0.0D0        
  368.       E = H            
  369.    10 IF (E .GE. Z) GO TO 20              
  370.          E = 2.0D0*E   
  371.          GO TO 10      
  372.    20 IF (E .LE. H) GO TO 30              
  373.          IF (E .LE. Z) Y = Y + E          
  374.          IF (E .LE. Z) Z = Z - E          
  375.          E = E/2.0D0   
  376.          GO TO 20      
  377.    30 Z = IDINT(Z + 0.5D0)                
  378.       Y = Y + Z        
  379.       IF (X .LT. 0.0D0) Y = -Y            
  380.       ROUND = Y        
  381.       RETURN           
  382.    40 ROUND = X        
  383.       RETURN           
  384.       END              
  385.       FUNCTION DFLOAT(I)
  386. C
  387. C   THIS IS THE AMIGA FUNCTION WHICH CONVERTS INTEGERS TO DOUBLE FLOATS
  388. C
  389.       IMPLICIT NONE
  390.       DFLOAT = DBLE(I)
  391.       RETURN
  392.       END
  393. SHAR_EOF
  394. #    End of shell archive
  395. exit 0
  396. -- 
  397. Bob Page, U of Lowell CS Dept.  page@swan.ulowell.edu  ulowell!page
  398. Have five nice days.
  399.